Tuning synchronization of integrate-and-fire oscillators through mobility 
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We analyze the emergence of synchronization in a population of moving integrate-and-fire oscil- 
lators. Oscillators, while moving on a plane, interact with their nearest neighbor upon firing time. 
We discover a non-monotonic dependence of the synchronization time on the velocity of the agents. 
Moreover, we find that mechanisms that drive synchronization are different for different dynamical 
regimes. We report the extreme situation where an interplay between the time-scales involved in 
the dynamical processes completely inhibit the achievement of a coherent state. We also provide 
estimators for the transitions between the different regimes. 
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After more than one decade of research on complex 
networks (HQ, where they have been understood as sim- 
ple collections of nodes and (weighted) links [1], time has 
come to analyze new paradigms of complex topologies. 
Examples of these new settings are, for instance, time 
dependent networks Q , networks at different layers (usu- 
ally called multiplex networks) [f|, spatial networks [f| 
and interdependent networks [7|. 

Concerning the case of time-dependent networks, up 
to now the interest has been focused on the need of find- 
ing new tools for characterizing topological properties of 
networks whose connectivity changes in time, but not as 
much on the interrelation between topological changes 
and dynamical evolution of the nodes themselves. This 
interrelation has been mainly considered for the asymp- 
totic cases in which characteristic times of both dy- 
namical processes are very different. When topological 
changes are very fast, agood proxy is the fast switching 
approximation (FSA) |8| whereas in the opposite limit, 
perturbation around the static case makes sense. In this 
paper we analyze dynamical effects in time dependent 
networks but understanding the network as the result of 
the motion of agents that interact when they are close 
enough[i{ni- 



When considering emergent properties of systems 
formed by elementary units with its own dynamics, one 
of the most important ones is synchronization [l2j ■ Actu- 
ally, synchronized behaviors appear in nature in groups 
that improve performance based on collaboration [13j 
and also in human actions it has shown to provide col- 
lective benefits [3, [l5| . The effect of changing patterns 
of interaction on synchronization features has been ana- 
lyzed in different settings, for instance in chemotaxis 11611 . 
mobile ad hoc networks [17| , wireless sensor networks [18[ 
and the expression of segmentation clock genes [l9| . 

Recently, a general framework of mobile oscillator net- 
works where agents perform random walks in a two- 
dimensional (2D) plane has been proposed [2C}. This 
framework, that reduces to FSA when velocity is high 
enough, is valid for models whose evolution can be well 
approximated by linear dynamics. This actually holds 
for models such as populations of Kuramoto oscillators 
[2l|,[22|, whose evolution, after a short transient time, is 



very well described by a set of linear equations that can 
be solved in terms of spectral properties of the Laplacian 
matrix [23|. Alternative approaches based on Fokker- 
Planck equations have been also proposed recently [24| . 

Here, we focus on a dynamical system, a population of 
integrate and fire oscillators (IFO), where the evolution 
takes place in two different time scales. One for the slow 
evolution of the internal state variables (the phase and 
the orientation) and the other for the instantaneous in- 
teraction between the units (pulse coupling). During the 
last years it has been shown that the interaction struc- 
ture plays a fundamental role in the dynamics of IFO 
networks [25[. Usually, IFO have been used to model 
neural systems but we can also find some examples of 
applications in other fields, as for example in economy 

In the present letter we consider a population of 
integrate-and-fire oscillators, which interact with their 
nearest neighbor only while freely displacing on a plane. 
Such a minimal interaction rule, in contrast to other ap- 
proaches based on an interaction radius 27 1. provides 
a strong non-monotonic dependence of the synchroniza- 
tion time on the velocity of the agents. Three different 
regimes arc identified according to the synchronization 
properties of the system: a slow regime, a fast switch- 
ing limit and an anomalous intermediate region between 
them. Remarkably, in this last region a divergence of 
the synchronization time is observed. Depending on the 
studied phenomena, synchronization can be taken either 
as as a positive or negative feature of the collective be- 
havior of the system. Global brain synchronization, for 
instance, is associated to epileptic seizures whereas local 
synchronization is related to some cognitive tasks [28j |. 
Synchronization in trading activities can be harmful for 
the stability of financial systems [l5[. Thus designing 
mechanisms that prevent global synchronization can be- 
come a very important tool in complex dynamical sys- 
tems 0. 

Our model consists on a population of N moving oscil- 
lators with identical velocity V and random orientation 
on a square of side length L with periodic boundary con- 
ditions. The internal phases of the agents (j> £ (0, 1) 
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increase uniformly with period r, 



d4n 
dt 



(1) 



until they reach a maximum value of 1, when a firing 
event occurs and the phase is reset. Upon such an event 
at time t, the firing oscillator influences its nearest neigh- 
bor (oscillator at minimal distance, labelled nn, see Fig. 
[T]4_) producing a random reorientation of its motion and 
an update in its phase by a factor e: 
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The phase update is performed at frozen time until the 
phases of all oscillators have been updated (some agents 
may reach their threshold and fire upon receiving a phase 
update from a firing neighbor). Then the phases evolve 
again uniformly in time until another update is triggered. 

The system is synchronized when a succession of con- 
secutive firing events (avalanche) equal to the system size 
N is detected. For the sake of clarity we define the (dis- 
crete) time T, as the number of times a given oscillator 
(that we identify with oscillator in our computer sim- 
ulations) has fired. This allows us to define T sync as 
the number of cycles this reference oscillator takes to en- 
ter the synchronized state (i.e., the number of updates 
needed for an avalanche of size N to occur). 

The chosen minimal interaction rule is such that the 
system lies far below the static percolation transition [3(| 
l3lT |. Therefore, global synchronization is not achievable 
without motion (since it is very unlikely that a giant 
connected component of size N exists HH ) . It is the non- 
null velocity of the oscillators that enables the system to 
reach the coherent state, therefore we could expect T sync 
(the average time the system needs to synchronize) to be 



a decreasing function of V, such that T s ., 



oo when 



V — > and T sync — > Tf > (a constant value) when V 
is high enough [35| , 

Fig. HP shows that, for V < V s , the synchronization 
time decreases as a power of V when the latter increases. 
Then, the decreasing slows down and T sync has a mini- 
mum at V = V m > V s . Beyond this value, the synchro- 
nization time gets larger and larger, until the system en- 
ters a region where it is unable to reach the coherent state 
in a finite time (gray area in Fig.[]Tj). For even larger val- 
ues of the velocity, V > V* , it decreases abruptly, finally 
reaching its asymptotic value when V = Vf. 

Three main regions are thus identified: A 'no- 
synchronization zone ' in the middle that separates a 'left 
region ' (small V) and a 'right region ' (high V). In the 
left region, we can separate two sub-regions: on the left, 
at V < V s , there is what we call the slow regime, and, 
on the right, a transition zone. The same happens for 
the right region. On the left, we find a transition region, 
while on the right, at V > V/, the system enters the fast 
limit, where T sync no longer depends on V. 
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FIG. 1: A Snapshot of the system in the incoherent state. 
Each disk is an oscillator, gray dashed arrows stand for first- 
neighbor relationships, while the continuous black one rep- 
resents a firing event that is taking place at that precise 
instant. B Final (T = T sync ) total interaction network at 
V = V m - Each node is an oscillator. Links represent shots, 
their weights are proportional to occurrence of the interac- 
tions among neighbors. Node color changes from purple to 
orange increasing the in-degree. Size increases with increas- 
ing out-degree. C The average synchronization time T 3ync as 
a function of V, for L = 400, N = 20, e = 0.1. Values of 
Vf and V s are calculated using estimators Q and (J6]l respec- 
tively. In the following, when not otherwise states, the values 
of the parameters are those used in this figure. Averages are 
performed over 2000 realizations and errorbars correspond to 
one standard derivation. 



Before entering into a more specific discussion we 
present a general characterization of the system behavior. 

In Fig. [31 upper panels, we introduce the cumulative 
individual interaction network (CIN) of an oscillator (la- 
belled 0), for two independent synchronization process at 
two different velocities (panel C: V = 2Vf while panel D: 
V = V m ). A CIN represents in a visual way the role 
played by a given unit in the signal spreading and is con- 
structed in the following way: whenever an oscillator i 
fires at oscillator and oscillator fires at oscillator j, a 
link between i and j is added. If the link already exists, 
its weight is increased. In the case does not receive any 
shot, we put a link between and j. A reciprocal shot 
is represented as a self-link. We repeat this process until 
the system reaches the synchronized state. 

To better understand the synchronization mechanisms 
we also report different relevant magnitudes. In Fig. [21 
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lower panels (panel A: V = 2Vf, panel B: V = V m ) we 
plot the global order parameter r](T) = (cos (27r<fo(T))) 
and in order to have an insight about what happens at lo- 
cal scale, the quantity A(T) = cos (27r0° n (T)) where </>° n 
is the phase of the unit to which the reference oscillator 
fired at time T. These parameters measure the synchrony 
of our system, ranging from a uniform phase distribution 
of our oscillators (r/(T) « A(T) « 0) to complete syn- 
chronization (rj(T) = A(T) = 1). Note that the averages 
are calculated with respect to all oscillators' phases {&} 
in the case of the global order parameter and only with 
respect to the nearest neighbor's phase <^ nn in the case 
of the local one [3(| . 

Finally, m(T) represents the total fraction of oscillators 
that have been out-neighbors of the reference oscillator 
up to time T. 

We observe strong differences between the two scenar- 
ios. When agents move fast, all nodes appear represented 
in the CIN since they play a global role sending and re- 
ceiving signals throughout the whole system. Interac- 
tions arc completely rewired at each time T, therefore 
m(T) increases very rapidly and <p nn is just a random 
variable extracted among N — 1 possible ones. This 
means that A is exactly the same as r], but with less 
statistics. Both quantities increase together (more or less 
noisily) because by means of the firing events the whole 
phase distribution becomes narrower. 

In contrast, when the mobility of the agents is reduced, 
few nodes conform the CIN. Each oscillator plays a local 
role mediating the interactions among a small number of 
units that are the same all the time, no matter how long 
the synchronization process could be. The behaviors of 
•q and A appear to be uncorrelated, being A = 1 almost 
all the time: each oscillator spends a long time with its 
neighbors, usually being able to synchronize with it be- 
fore changing. At the beginning (T < 500) whenever it 
starts firing towards a new oscillator (black vertical lines 
in Fig. [2]) A experiences an abrupt decreasing while m(T) 
increases as it is the first time the reference oscillator 
meets that neighbor. Later (T > 500), the chances to 
change a neighbor for another one already known, hav- 
ing a very similar phase, increase. The phase distribution 
becomes narrower and specially at local scale, among the 
units the oscillator of reference can meet, the dispersion 
is small. Consequently neighbor changes do not affect A 
anymore. 

The observations point out that the fast regime can 
be understood as an homogeneous regime while slow ve- 
locity enhance heterogeneity among units. Hence the 
mechanisms that allow the system to synchronize have 
to be different in the two cases. The system has different 
strategies to reach the coherent state in the left and in the 
right region but neither of them work in the intermediate 
region. 

In the following paragraphs we will determine the 
region of the parameters space corresponding to each 
regime providing a quantitative estimator for both V s 
and Vf as functions of N, e, and L. 
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FIG. 2: Upper panel: Final (T = T sync ) network of the 
interactions mediated by a single oscillator (labelled "0"), in 
the fast limit, at V = 2Vj (panel A) and at V = V m (panel 
B) respectively. Node and link color and size code as in fig- 
ure 1 B. Lower panels: Parameters n (blue), A (red) and 
m (black) are plotted as a function of time, from T = to 
the synchronization time for a given oscillator and a single 
realization of the process. Yellow vertical lines mark a change 
of the in-neighbor, while the purple lines stand for a change 
of the out-neighbors. In panel A, the velocity is V = 2V/. 
In panel B, it is V = V m - Where yellow lines are so dense 
that they form a yellow band, it means that the oscillator has 
more than one in-neighbor simultaneously. 



In order to understand the different regimes we need 
to compare the time scales of the dynamical processes 
that lead the system to its final synchronized state, i.e. 
how fast arc oscillators to change neighbors effectively 
and how fast they are able to locally synchronize. 

Let us first focus on the fast regime. As a kind of 
mean-field assumption we can imagine a single oscilla- 
tor, i, moving on a plane where the rest of the popu- 
lation is fixed in the most disperse possible configura- 
tion, a square lattice such that the system is divided into 
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(N — 1) squares of side length £ = L/y/N — 1. Oscillator 
i changes its neighbor when it exits a square to enter one 
of the adjacent ones. It will do so between time T — 1 and 
T if, at time T — 1, the component of V perpendicular to 
that side is larger than the distance s separating i from 
that same boundary. First, notice that if V > £ then the 
exit probability is equal to 1 and each unit changes its 
neighbor at each time step. One can conclude that £ is a 
good estimator for the value of the velocity at which the 
system enters the fast switching limit, i.e. 



(3) 



On the opposite limit, when oscillators move very 
slowly (VVN /L <C 1), we can also compute the average 
time spent in changing neighbors (see the Supplemental 
Material [32|) 



Wy/N - 1 



(4) 



This time needs to be compared with the time a char- 
acteristic cluster needs to locally synchronize [2(j ■ Given 
the minimal rule of interaction we have proposed, the 
average size of (weakly) connected components in the re- 
sulting network of our model + 3 v / 3)/3 
(see [1U and [IH). For such a characteristic cluster of 
size 3 there is only one possible configuration [37| (see 
Fig. [T]4.) whose synchronization time is given empirically 
by the expression 



T 3 = n/e. 



(5) 



with k = 2.0 ±0.1 [32j. A necessary condition for the sys- 
tem to reach complete local synchrony is that, on average, 
no topological change has to occur during a time frame 
T3. Since the average time separation between two neigh- 
bor changes of the same unit is T out , the average time 
separation between two topological changes in the system 



can be expressed as T ( 
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T out /N~TrL/(WN 3 / 2 ). 



With these two ingredients, we can easily identify the 
starting point of the slow regime as that set of parameters 
values such that T 3 
we obtain 



T^. ■ Hence, fixing N, e and L, 



V, 



7r eL 
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Both Vf and V s are reported in Fig.Q] and have been 
checked for various values of e, N, and L, showing a very 
good agreement. In particular, in Fig.[3]we show how our 
predictor mark the transition point of the slow regime. 
Moreover, we can identify a universal minimum for all 
the curves marking the value of V m . 

Summarizing, we have proposed a model of moving 
integrate-and-fire oscillators that interact only with the 
nearest neighbor. This minimal interaction rule can ac- 
count for some physical situations where communication 
between agents is minimized. We have found three very 
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FIG. 3: Average T sync as a function of the ratio Tz/T^^ for 
several values of TV (filled squares: TV = 10, filled circles: TV = 
20, triangles: TV = 30, empty squares: TV = 40, empty circles: 
TV = 50) and e. The considered values of e are, from the top to 
the bottom: e = 0.05 (orange), e = 0.1 (red), e = 0.2 (black); 
additionally, for TV = 20, we also considered e = 0.02 (yellow 
top curve) and e = 0.3 (purple bottom curve). The vertical 
dashed line at Ts/T^ = 1 corresponds to V = V s , the other 
one to V — Vm- Averages have been performed over 2000 (for 
N < 20) or 1000 (N > 20) realizations. 



different dynamical regimes. In the region where oscilla- 
tors move very slowly and oscillators tend to keep the 
same neighbors for a long time, local synchronization 
dominates and a global one proceeds through a very large 
number of interactions (firings and changes of neighbors). 
In the other limit, of very fast motion, oscillators fire 
at random neighbors which makes that interactions to 
be very effective and phases approach each other in a 
monotonous way and synchrony is achieved directly on 
the global scale. Interestingly, we find an intermediate 
regime where none of these mechanisms work and the sys- 
tem is not able to reach synchronization. The presence 
of a no synchronization band needs to be further studied 
and its robustness checked for other types of non-linear 
interacting oscillators in future works. 

Despite synchronization is usually seen as a positive 
outcome of a cooperative dynamical system, under some 
circumstances, as it happens, for instance, in brain dy- 
namics, global synchronization is not a healthy state. 
Our model represents a paradigmatic example of a sys- 
tem where synchronization can be prevented by tuning 
the rate at which the topology changes. We have deter- 
mined the bounds of this region in terms of the charac- 
teristic times of the model, those related to local syn- 
chronization and to mobility. 
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